home *** CD-ROM | disk | FTP | other *** search
/ Compendium Deluxe 2 / LSD and 17bit Compendium Deluxe - Volume II.iso / a / prog / asmsrc / thesource-7.lha / Source / polyfit.lha / polyfit / mat.c < prev    next >
C/C++ Source or Header  |  1990-05-04  |  1KB  |  76 lines
  1. #include "fit.h"
  2.  
  3. matrix MatInit(n)
  4.    register int    n;
  5. {
  6.    register int    i;
  7.    register matrix mat;
  8.  
  9.    mat = (matrix)Calloc(n, sizeof(matrow));
  10.  
  11.    for (i = 0;  i < n;  i++)
  12.       mat[i] = (matrow)Calloc(n + 1, sizeof(matelm));
  13.  
  14.    return (mat);
  15. }
  16.  
  17.  
  18. double MatInv(n, mat)
  19.    register int    n;
  20.    register matrix mat;
  21. {
  22.    int    pivRow;
  23.    double piv, det, div, pivTmp;
  24.    matrow tmpRow;
  25.    register int    p, r, c;
  26.  
  27.    for (det = 1.0, p = 0;  p < n;  p++)
  28.    {
  29.       for (piv = mat[pivRow = p][p], r = p;  r < n;  r++)
  30.          if (((pivTmp = mat[r][p]) < 0 ? -pivTmp : pivTmp) > piv)
  31.             piv = mat[pivRow = r][p];
  32.  
  33. #ifndef NO_DETERM
  34.       if ((det *= piv) == 0.0)
  35.          break;
  36. #endif
  37.  
  38.       if (pivRow != p)
  39.       {
  40. #ifndef NO_DETERM
  41.          det = - det;
  42. #endif
  43.          tmpRow = mat[p];
  44.          mat[p] = mat[pivRow];
  45.          mat[pivRow] = tmpRow;
  46.       }
  47.  
  48.       for (c = p;  c <= n;  c++)
  49.          mat[p][c] /= piv;
  50.  
  51.       for (r = 0;  r < n;  r++)
  52.          if (r != p)
  53.             for (c = n;  c >= p;  c--)
  54.                mat[r][c] -= mat[r][p] * mat[p][c];
  55.    }
  56.  
  57. #ifndef NO_DETERM
  58.    return (det);
  59. #endif
  60. }
  61.  
  62.  
  63. pointer Calloc(n, s)
  64.    unsigned n, s;
  65. {
  66.    register pointer ret;
  67.  
  68.    if ((ret = (pointer)calloc(n, s)) == NULL)
  69.    {
  70.       perror("Allocating matrix");
  71.       exit(1);
  72.    }
  73.  
  74.    return (ret);
  75. }
  76.